#ifndef  __PARAMETERS_H_
#define  __PARAMETERS_H_ 

#include <cmath> // sqrt
#include <iostream>
#include <map>
#include <vector>

#include "tppmc/Enumerations.h" // For LHAPDF::Flavour

using namespace std;

#include "tppmc/Particle.h"

// Partonic content - for transversity code 
struct PARTONCONTENT {
   double
   up,
   down,
   anti_up,
   anti_down,
   strange,
   anti_strange,
   charm,
   anti_charm,
   bottom,
   anti_bottom,
   top,
   anti_top,
   glu;
   
};

/**
 This is a temporary treatment to make it easier to bridge old and
 new implementation.
 Needs to be reworked.
 */
class FitValues {
   
public:
   
   FitValues(double a = 0., double b = 0., double n = 0., double m2 = 0.)
   : mA(a), mB(b), mN(n), mM2(m2) { }
   
   virtual ~FitValues() { }
   
   double Alpha() const { return mA; }
   double Beta() const { return mB; }
   double N() const { return mN; }
   double M2() const { return mM2; }
   double M() const { return ::sqrt(M2()); }
   
   void Speak() const {
      std::cout << mA << "\t" << mB << "\t" << mN << "\t" << mM2 << std::endl;
   }
   /**
    Evaluate the result of this fit at the x of a parton.
    */
   double EvaluateTransversity(const tppmc::pythpart&) const;
   /**
    Evaluate the result of this fit at the x and kT of a parton.
    */
   double EvaluateSivers(const tppmc::pythpart&) const;
   
protected:
   
   double mA;
   double mB;
   double mN;
   double mM2;
};

// Parameters class
class VALUES {
private:
public:
   VALUES() {};
   ~VALUES() {};
   
   // Down quark
   double a_down;
   double b_down;
   double n_down; 
   double m2_down; // kt dependence
   
   
   // Up quark
   double a_up;
   double b_up;
   double n_up;
   double m2_up; // kt dependence
   
   // anti Down quark
   double a_anti_down;
   double b_anti_down;
   double n_anti_down; 
   double m2_anti_down; // kt dependence
   
   
   // anti Up quark
   double a_anti_up;
   double b_anti_up;
   double n_anti_up;
   double m2_anti_up; // kt dependence
   
   // Strange quark
   double a_strange;
   double b_strange;
   double n_strange; 
   double m2_strange; 
   
   
   // anti Strange quark
   double a_anti_strange;
   double b_anti_strange;
   double n_anti_strange;
   double m2_anti_strange;
   
   // Quick work-around to give subscripted access to parameters
   // based on LHAPDF::Flavour.
   // Gradually phase out the above individual double values in
   // favour of this.
   
   std::map<int, FitValues> mValues;
   
   FitValues* Values(int);
   const FitValues* Values(int) const;
   
   // Temporary...
   void Populate();
};

const int set_number = 200; // we will generate 200 sets

class SETS{ // organize all sets here params, errors, sets of params
private:
   VALUES main_parameters;
   const char* FILE_READ; // file to read parameters
   const char* SETS_FILE_READ; // file to read set parameters
   const char* FILE_WRITE; // file to write parameters
   const char* SETS_FILE_WRITE; // file to read set parameters
   
public:
   SETS() {}; // constructor...
   ~SETS() {}; // deconstructor...
   
   VALUES parameters;
   VALUES errors;
   
   VALUES set[set_number];
   void SetReadFile(const char* file);
   void SetSetsReadFile(const char* file);
   void SetWriteFile(const char* file);
   void SetSetsWriteFile(const char* file);
   void ReadParameters();
   void ReadSets();
   void WriteParameters();
   void WriteSets();
   void UseMainSet(void);
   void UseSet(int k);
   
   vector<double> MinMax( double (*f)(double x), double x ); // evaluate min and max of a function f(x) on the parameters sets
};

// HERE WE DEFINE ALL PARAMETERS...
class PARAMETERS{
private:
   // All params...
   double kt2_average;
   double ptq2_average;

public:
   
   PARAMETERS(); // constructor...
                 //Get params...
   double GetKt2Average( void ) const { return kt2_average;};
   double GetPtq2Average( void ) const { return ptq2_average;};

   /**
    Print all parameters to an output stream.
    */
   virtual void Speak(std::ostream&) const;
   
   SETS Sivers;
   SETS Transversity;
   SETS Collins;
};

#endif // #ifndef  __PARAMETERS_H_
